#!/usr/bin/env python3
"""
Calculate Pi Using Various Methods

This script implements several algorithms for estimating the value of π (pi):
1. Monte Carlo Method - using random sampling
2. Leibniz Formula - using infinite series
3. Bailey-Borwein-Plouffe Formula - allows calculating specific digits
4. Nilakantha Series - converges faster than Leibniz
5. Archimedes' Method - using polygons inscribed in a circle

Each method demonstrates a different mathematical approach to approximating π.
"""

import random
import math
import time


def monte_carlo_pi(num_samples=1000000):
    """
    Estimate π using the Monte Carlo method.
    
    This method uses random sampling to estimate π by comparing the ratio of 
    points falling inside a circle to those falling inside a square.
    
    Args:
        num_samples: Number of random points to generate
        
    Returns:
        Estimation of π
    """
    inside_circle = 0
    
    for _ in range(num_samples):
        # Generate random point in a 1x1 square
        x = random.random()
        y = random.random()
        
        # Check if point is inside the circle
        if x**2 + y**2 <= 1:
            inside_circle += 1
    
    # Ratio of area of circle to square is pi/4
    return 4 * inside_circle / num_samples


def leibniz_pi(iterations=1000000):
    """
    Calculate π using the Leibniz formula.
    
    Uses the infinite series: π/4 = 1 - 1/3 + 1/5 - 1/7 + 1/9 - ...
    
    Args:
        iterations: Number of terms in the series to sum
        
    Returns:
        Estimation of π
    """
    pi = 0.0
    sign = 1
    
    for i in range(iterations):
        pi += sign * (1.0 / (2 * i + 1))
        sign *= -1  # Alternate sign
    
    return 4 * pi


def bbp_pi(n_digits=100):
    """
    Compute π using the Bailey-Borwein-Plouffe formula.
    
    This formula allows calculating specific digits of π without
    calculating previous digits.
    
    Args:
        n_digits: Number of iterations for the computation
        
    Returns:
        Estimation of π
    """
    pi = 0
    for k in range(n_digits):
        pi += (1.0/16**k) * (
            4.0/(8*k+1) - 
            2.0/(8*k+4) - 
            1.0/(8*k+5) - 
            1.0/(8*k+6)
        )
    return pi


def nilakantha_pi(iterations=10000):
    """
    Calculate π using the Nilakantha series.
    
    Uses the series: π = 3 + 4/(2×3×4) - 4/(4×5×6) + 4/(6×7×8) - ...
    This series converges faster than the Leibniz formula.
    
    Args:
        iterations: Number of terms in the series to sum
        
    Returns:
        Estimation of π
    """
    pi = 3.0
    sign = 1
    
    for i in range(1, iterations + 1):
        denominator = (2 * i) * (2 * i + 1) * (2 * i + 2)
        pi += sign * (4.0 / denominator)
        sign *= -1  # Alternate sign
    
    return pi


def archimedes_pi(n_sides=1000000):
    """
    Estimate π using Archimedes' method with n-sided polygons.
    
    This geometric approach uses polygons inscribed in a circle
    to approximate π.
    
    Args:
        n_sides: Final number of sides in the polygon
        
    Returns:
        Estimation of π
    """
    # Start with a hexagon (n=6)
    perimeter = 6.0  # For a unit circle
    current_sides = 6
    
    while current_sides < n_sides:
        # Double the number of sides
        current_sides *= 2
        
        # Calculate new perimeter based on geometric properties
        # This formula calculates the side length of a polygon with 
        # twice as many sides, given the perimeter of the current polygon
        side_length = math.sqrt(2 - 2 * math.sqrt(1 - (perimeter / (2 * current_sides / 2))**2))
        perimeter = current_sides * side_length
    
    # π ≈ perimeter / 2 (for a circle of radius 1)
    return perimeter / 2


def main():
    """Test and compare different methods for calculating π."""
    true_pi = math.pi
    
    print(f"True value of π: {true_pi}\n")
    print("Testing different methods to estimate π:\n")
    
    # Dictionary to store results for comparison
    results = {}
    
    # Monte Carlo method
    start_time = time.time()
    mc_pi = monte_carlo_pi(1000000)
    mc_time = time.time() - start_time
    results["Monte Carlo"] = (mc_pi, mc_time)
    
    # Leibniz formula
    start_time = time.time()
    leibniz_result = leibniz_pi(1000000)
    leibniz_time = time.time() - start_time
    results["Leibniz"] = (leibniz_result, leibniz_time)
    
    # BBP formula
    start_time = time.time()
    bbp_result = bbp_pi(10)
    bbp_time = time.time() - start_time
    results["BBP"] = (bbp_result, bbp_time)
    
    # Nilakantha series
    start_time = time.time()
    nilakantha_result = nilakantha_pi(10000)
    nilakantha_time = time.time() - start_time
    results["Nilakantha"] = (nilakantha_result, nilakantha_time)
    
    # Archimedes' method
    start_time = time.time()
    archimedes_result = archimedes_pi(96)
    archimedes_time = time.time() - start_time
    results["Archimedes"] = (archimedes_result, archimedes_time)
    
    # Print results in a formatted table
    print(f"{'Method':<15} {'Estimated π':<20} {'Error':<15} {'Time (s)':<10}")
    print("-" * 60)
    
    for method, (value, exec_time) in results.items():
        error = abs(value - true_pi)
        print(f"{method:<15} {value:<20.15f} {error:<15.10e} {exec_time:<10.6f}")


if __name__ == "__main__":
    main()